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Abstract 

During the past five years the low temperature heat capacity of simple semiconductors and 
insulators has received renewed attention. Of particular interest has been its dependence on isotopic 
masses and the effect of spin- orbit coupling in ab initio calculations. Here we concentrate on the 
lead chalcogenides PbS, PbSe and PbTe. These materials, with rock salt structure, have different 
natural isotopes for both cations and anions, a fact that allows a systematic experimental and 
theoretical study of isotopic effects e.g. on the specific heat. Also, the large spin-orbit splitting of 
the 6p electrons of Pb and the 5p of Te allows, using a computer code which includes spin-orbit 
interaction, an investigation of the effect of this interaction on the phonon dispersion relations 
and the temperature dependence of the specific heat and on the lattice parameter. It is shown 
that agreement between measurements and calculations significantly improves when spin-orbit 
interaction is included. 

PACS numbers: 63.20.D-, 65.40.Ba 
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I. INTRODUCTION 



Considerable effort has been spent recently in the investigation of the heat capacity 
C of semiconductors and insulators in the region below the Debye temperature Or, , ui 
particular around 0.1 6'/) , where strong deviations from Debye's power law take place. 
The availability of stable isotopes has also enabled the investigation of the dependence of C 
on one isotopic mass (in monatomic crystals such as C (Ref. Q), Si (Ref. Ge (Ref. 130), Sb 
(Ref. |5|), Bi (Ref. |6|y) or on the isotopic mass of each constituent in polyatomic compounds, 
such as GaN (Ref. |7|) or ZnO (Ref. |8|), whereas the development of efficient computer codes 
for electronic band structure calculations has made possible ah initio calculations of the 
phonon dispersion relations and the temperature dependence of the specific heat. 

Although experimental and calculated results for ~ Cp (Ref. jo]) agree usually reason- 
able well, it has been recently discovered that considerable discrepancies exist when heavy 
constituent atoms are present, e.g. Bi (Ref. Pb in PbS (Ref. , if the Hamiltonian used 
for the calculation of the electronic structure does not include spin-orbit (s-o) interaction. 
By performing ab initio calculations in which the s-o interaction is switched on and off, it 
has been rather conclusively shown that this interaction "softens" the phonon frequencies 
(Ref. and thus increases the low-temperature maximum found in Ci,/T^ versus T.— 

In this article we present ab initio calculations of the phonon dispersion relations of 
PbS, PbSe and PbTe based on the electronic band structure obtained with the ABINIT 
code.— This program enables the inclusion of s-o interaction, separately or jointly for the 
cation or anion constituents, thus making it possible to separate the s-o contributions to the 
dispersion relation and the heat capacity. The ABlNlT program determines the dynamical 
matrix elements by perturbation theory. Once these elements are known, it is a simple 
task to calculate the dispersion relations for a given set of isotopic masses and thus to 
determine isotopic effects on the dispersion relations and the specific heat. Because of the 
importance of s-o interaction, we have performed the isotope effect calculations only with s-o 
interaction included. Within the range of stable isotopes available in nature (and accessible 
to our budget) only mass changes of a few percent are possible: the resulting variations in 
the physical properties are thus linear in the mass changes. Results for PbS have already 
been published in Refjoj. At that time, the ABINIT code did not properly include s-o 
interaction for diatomic polar compounds. We therefore present here similar results including 
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s-o interaction: the inclusion of this interaction reduces the discrepancy between calculated 
and experimental values of C^^p (29%) by a factor of 2. A similar reduction is found for PbSe, 
whereas for PbTe the calculated values almost exactly agree with the measured data. In 
this latter case the effect of the s-o interaction on the maximum of C^,/T^ is 18% (we recall 
that the s-o interaction softens the phonons and correspondingly increases C^/T'^ at the 
maximum). Globally, the inclusion of s-o interaction also improves the agreement between 
the calculated and the measured (by inelastic neutron scattering, INS) phonon dispersion 
relation, except for the LO phonons in the vicinity of the F-point of the Brillouin zone (BZ). 
At this point, the long range electric fields associated with the LO phonons lead to some 

lich we have not been able to avoid completely. 



convergence problems w 

In Ref. y, and in Ref . |13| we also reported densities of one- and two-phonon (with the same 
/c-value both, as required to compare with optical spectroscopy results) states, calculated for 
PbS in the absence of s-o interaction. We present here similar two-phonon spectra obtained 
for the three lead chalcogenides including s-o interaction. The features obtained in the 
density of one-phonon states, and their projections on the constituent atoms, are useful 
(and very instructive) for the interpretation of the isotope effects. 

We have also looked at the effect of s-o interaction on the lattice parameter as obtained 
by energy minimization. This effect is rather small (~ 0.2%). It does not help to reduce the 
discrepancy between the calculated and the measured values, which is about 2%. 

II. THEORETICAL DETAILS 

The calculations of the dispersion relations and the specific heat were performed with the 
dynamical matrix obtained from the LDA electronic structure calculations using the ABINIT 
code.— Hartwigsen-Goedecker-Hutter relativistic separable dual-space pseudopotentials were 
used.— We performed checks with LDA and GGA functionals and convinced ourselves that 
no significant differences resulted; we therefore used LDA functionals throughout. For the 
Pb pseudopotential we also checked that our conclusions did not depend on the inclusion of 
bd electrons in the valence bands. This allows us to use the pseudopotential as it is given 
in the ABINIT website. The implementation of the s-o term has been discussed in Ref. lU 



We investigated cell parameter convergence as a function of energy cutoff and /c-grid mesh. 
With the cutoff at 60 Hartree the total energy is converged up to 0.5 meV and stresses 
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are lower than 0.006 GPa. The BZ is sampled using a 6 x 6 x 6 Monkhorst-Pack grid.— 
Technical details of the calculations of the phonon dispersion relations can be found in Ref. 
and 
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171 . Prior to obtaining the dynamical matrix the lattice parameter was optimized 
through minimization of the total energy. In the calculations we used this parameter and 
not the experimental one. 



After obtaining the dispersion re^ 
the integral given in Eq. (2) of Ref 



.y. 



ations, the phonon free energy F was calculated with 
The specific heat was obtained with the expression: 

c.^-r(^)^. (1) 

After calculating the dynamical matrix elements we diagonalized the Hamiltonian for 
two different sets of masses (either cation or anion) differing by about 5%. The resulting 
dispersion relations were then used to calculate the derivatives of C„ with respect to either 
mass. The logarithmic derivatives were also calculated because they can be related to the 
corresponding derivative with respect to temperature. The appropriate relation is given in 
Eq. (4) of Ref. O for the case of a monatomic crystal. For diatomic crystals such relation 
is not as simple since both, derivatives of versus the masses of cation and anion must be 
added in order to obtain the temperature derivative: 

rflnMpb dlnMx 2 ^ cilnT ^ ' 

where X = S, Se or Te, respectively. 

When calculating numerically the derivatives above at very low temperatures (< 5K) one 

encounters a convergence problem because of the small values of and T. Fortunately, 

the limit T ^ can be obtained analytically by usingthe Debye-T"^ approximation. The 

corresponding expressions are found in Eq. (5) of Ref. [9. 



III. EXPERIMENTAL PROCEDURE 



Samples of PbX (X = S, Se, Te) were prepared by first reacting the corresponding pure 
elements and then subliming the product in an argon atmosphere. 

In order to purify the lead isotope oxide layers on the metal pieces were removed by 
etching in diluted nitric acid and, in the case of natural lead, by cutting. Subsequently, the 
pieces were melted in silica ampoules in hydrogen or argon. Remaining oxide stuck to the 
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ampoule wall after rocking and rolling the droplets at temperatures from 400 to 650°C. The 
chalcogens were purified by sublimation and separation of the ampoule portion containing 
the residues. In the case of ^^'^Te, we first reduced ^^°TeO with sulfur.— The synthesis was 
performed in argon by increasing the temperature to ~650 °C of the chalcogen. In the 
sealed silica ampoule crystal growth took place by sublimation at temperatures of 750 to 
850 °C under excess of the corresponding chalcogen. This resulted in a cabbage-like growth 
of the compound on the lead. A few platelets up to 4 x 4 mm^ were obtained during 1 to 
2 weeks.— The preparation and some properties of the PbS samples are described in some 
of the measured PbS samples were natural galena crystals.— 

The heat capacities were measured between 2 and 280K with a PPMS system (Quantum 
Design, 6325 Lusk Boulevard, San Diego, CA.) as described in detail Ref. |8. 

We also measured natural galena crystals. 



IV. DISPERSION RELATIONS 



The phonon dispersion relations calculated along three high symmetry directions of the 
BZ ( [111], [100], [Oil]) for PbS are shown in Fig. □ together with INS data of Elcombe.-^i 
The calculations without s-o coupling were already published in Ref. O. Note that the s-o 
interaction lowers all phonon frequencies, bringing the calculations in better agreement with 
the experiments except for the LO band where, on the average, the agreement is similar. As 
we shall see later, the decrease in the calculated frequencies raises he maximum in C^/T^ , 
thus improving agreement with experimental data. 

Figure [2] and [3] display the dispersion relations of PbSe and PbTe respectively, calculated 
with and without s-o coupling, together with experimental data. For both materials we 
also observe the decrease of phonon frequencies when the s-o interaction is taken into ac- 
count. For PbSe only few experimental (INS) data points are available.— Figure [2] provides 
some indication that also in this case the calculations with s-o coupling come closer to the 
experimental data than those without. 

Figure [3] shows calculated (with and without s-o interaction) and measured (by INS at 
297K) phonon dispersion relations of PbTe. Again we observe a lowering of nearly all fre- 
quencies induced by taking into account s-o interaction, which also brings the calculations, 
in most cases, closer to the experimental data. The measured TO frequencies are known to 
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FIG. 1: Phonon dispersion relations of PbS with natural isotope composition of Pb and S calculated 
with and without s-o split ting within the harmonic approximation. The points were obtained by 



INS at 300 K (From Ref. 
-1 



21 



). Typically, anharmonic effects should have lowered them by ~ 2 



cm 



be strongly renormalized downwards by the anharmonic interaction, even at low tempera- 



tures. By linearly extrapo 
temperature (from Ref. 



, ating to T=0 the measured dependence of the TO frequency on 
25l ) we obtain an unrenormalized (harmonic) TO frequency of ~13 
cm~^, rather close to that calculated with s-o interaction (~ 14 cm~^, cf. FigjS]). This low 
TO frequency for T ~ OK underscores the nearly ferroelectric nature of this material which 
would nominally become ferroelectric at ~-70K.— 

A disturbing feature of Fig [3] is the dip in the LO-phonon band around the F- point when 
s-o coupling is taken into account. We have been able to show that the appearance of this 
dip is related to the closing and sign reversal of the electronic gap induced by the s-o splitting 
of the valence and conduction band edges. Lowering the s-o Hamiltonian by multiplication 
by an adjustable parameter A < 1 (cf. Ref. ^, see also Sect. IV) the electronic gap of PbTe 
decreases from 0.32 eV for A = 1, becomes zero at A ~ 0.75 rising again to 0.53 eV for A = 
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FIG. 2: (color online) Phonon dispersion relations of PbSe with natural isotope composition of Pb 
and Se calculated with and without s-o splitting within the harmonic approximation. The circles 
(o) were obtained by INS at 300 K (From Ref. |22| ). The (red) square was obtained by tunnel 



spectroscopy at 4.2 K.— Typically, anharmonic effects should have lowered them by ~ 2 cm~^. 
The (red) closed circle was obtained by ir transmission at 1.4K.— The (red) dotted line isuggests 
how the LO band should look like (for more details see text). 



(the value of the gap measured at 4 K is 0.19 eV (Ref. l27|). Note that for PbSe (Fig.2) 
a similar but smaller dip occurs when the s-o interaction is switched on. In this case, the 
band gap nearly vanishes for A = and becomes 0.52 eV for A = 1, with the same inverted 
order as in PbTe. Reversal of the gap is smaller than for PbTe because the s-o splitting of 
Se is smaller than that of Te. For PbS (cf. Fig. [T]) no dip is present other than that already 
existed without s-o splitting. The gap is ~ for A = and 0.47 eV (also inverted) for A =1. 

In view of the systematics of the effects just described we conjecture that the anomalous 
dip in the LO frequencies observed at the F-point, is due to electron-phonon interaction of 
the type named after Kristoffel and Konsin.— This interaction lowers the phonon frequencies 
at F and can even lead to ferroelectricity when a phonon frequency vanishes (this is almost 
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FIG. 3: (color online) Phonon dispersion relations of PbTe with natural isotope composition of Pb 



and Te calculated with and without s-o splitting wit 
circles were obtained by INS at 297 K (From Ref 



lin the harmonic approximation. The (black) 



26l ). The (red) squares were obtained by optical 



spectroscopy at 300K (32 cm ^) and 5K (18 cm ^).— The (red) dotted line suggests how the LO 
band should look like (for more details see text) 

but not quite the case for PbTe but happens for SnTe, for which the gap is reversed with 
respect to PbTe.— Since the measured LO frequency does not support the presence of the 
dip calculated for PbTe (and also less strongly for PbSe) we conclude that the dip is due to 
some as yet unidentified effect of the s-o interaction, probably related to the gap anomalies 
just described (and to the so-called gap problem ubiquitous in LDA calculations). We have 
just made a guess as to how the LO band should look like once the problem is removed, and 
represented it by a dotted red line, for both PbTe and PbSe (cf. Fig. [2] and [3]). Because 
of the small volume of A;-space encompassed by the dip (~ 0.1% of the BZ) we believe that 
this problem should not affect our calculations of heat capacities. 
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V. PHONON DENSITIES OF STATES AND RAMAN SPECTRA 



The densities (DOS) of one-phonon states for the three materials under consideration are 
displayed in Fig. HI normalized to 6 states, as corresponds to one primitive cell (PC). This 
figure shows not only the total DOS but also its projection onto each of the two constituent 
atoms. A DOS calculated usingthe dispersion relations without s-o interaction has already 
been reported for PbS in Ref. |9|. Here we report only the calculations performed with s-o 
interaction. 

A comparison of the results in Fig. H] for PbS with those in Fig. [2] of Ref. 1^ (no s-o 
interaction) clearly shows the shrinking of the width of the phonon bands, from 250 cm~^ 
without s-o to ~ 235 cm"^ produced by the s-o interaction. Also, whereas without s-o 
interaction the vibrations of the two constituents are almost fully separated in frequency 
(dividing line 112 cm~^), this is not the case any more if s-o coupling is included. As shown 
in Fig. Hlfor PbS both types of vibrations (Pb-like and S-like) overlap in the region between 
100 and 130 cm~^, presumably related to the shrinking of the phonon band width. In the 
case of PbSe three bands are seen in Fig. HI The low frequency one, with a strong peak at 
~ 50 cm~^, corresponds to Pb-like vibrations of TA phonons, whereas the high frequency 
one, centered at about 125 cm~^, is LO- and TO- like and encompasses almost exclusively 
(95%) Se vibrations. The third band, centered around 80 cm~^, is a mixture of LA- and 
TO- like modes and encompasses vibrations of both constituent atoms with nearly the same 
amplitude. The DOS of PbTe also exhibits 3 bands with a "pseudogap" at 60 cm~^ and 
a second, not as deep, at ~ 85 cm~^ .The lower band, mainly TA-like with some TO and 
LA contribution, is predominantly Pb-like, whereas the middle band is about 66% Te-like 
and 33% Pb-like. These facts, which will be useful when we discuss the dependence of the 
heat capacity on isotopic mass, are quantitatively represented in Fig. |4] (right vignette) 
which displays the integrated total number of states and their projections, as a function of 
frequency. The integrated (total) phonon states of PbS show clearly the almost complete 
separation of Pb and S vibration, below 115 cm~^ Pb-like, above S-like, and the existence 
of a central band of mixed modes in PbSe and PbTe. 

Beside the one-phonon DOS just described, we have also calculated the DOS of two- 
phonon states with total zero wave vector (or, equivalently, equal wave vectors). These DOS 
are useful for the interpretation of second order Raman scattering since first order Raman 



10 




FIG. 4: (color online) (left) (a-c) Phonon density of states PDOS of PbS, PbSe and PbTe, from 
top to bottom, respectively. The projections of the PDOS on the two constituents is shown by 
the dashed lines, red: Pb projection; blue: projection on S, Se and Te, respectively, (right) (a-c) 
Integrated PDOS (solid black line) showing the total number of states and their projections on the 
cation (red line) and the anion (blue line) constituents. 
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FIG. 5: (color online) (a-c) Densities of two-phonon states (sum and difference) for PbS, PbSe, 
and PbTe, from top to bottom respectively, with zero total /c-vector. Notice that for the sums the 
curves are normalized to 36 states (6x6 states per PC). For the differences they are normalized 
to 15 states (6x5/2 staes per PC). The blue arrows mark structure which is weakly (except for 
strong 430 - 460 cm~^peaks) observed in the Raman spectrum of PbS [ at ~ 150, 180, 205, 220, 



250, 320, 340, 430, 460 cm 
published in Refs. 131 



and 
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The green arrows indicate possible structure related to the spectra 



32 for PbSe. 



scattering is forbidden, cf. Ref. Il3|. Two kinds of second order processes are possible, one in 
which the frequencies of the two phonons add (sum processes) and the other in which they 
subtract (difference processes). We calculate the DOS for both. In Raman scattering sum 
processes are present even at T=0 and increase with increasing temperature according to 
the appropriate statistical factors. Difference processes have vanishing intensity at T=0 and 
can only be observed with increasing temperature. We display in Fig. [5] the densities of two- 
phonon states for such sum and difference processes. Note that they are normalized to 36 
states (6x6 states per PC) for the former and to 15 states (6x5/2) for the latter, the factor 
of 2 corresponding to the separation of Stokes and anti-Stokes spectra. The blue arrows 
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mark structure which is weakly observed (except for the strong 430 - 460 cm ^ peaks) in 
the Raman spectrum of PbS [ at ~ 150, 180, 200, 220, 250, 320, 340, 430, 460 cm'^.i^ The 
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and 



32 



green arrows indicate possible structure related to the spectra published in Refs. 
for PbSe. No structure has been identified as corresponding to phonon differences, although 
the band seen in Fig. [Slat ~ 40 cm~^ , which seems to correspond to differences of TO and 
LA phonons (see Fig. [1]), should increase strongly with temperature and be identifiable with 
a good spectrometer that covers, with little straight light, the corresponding region. 

Second order Raman spectra of PbSe and PbTe are rather scarce in the literature. Some 



data are available for PbSe films grown on BaF2 (Ref. l3ll ) and for bulk PbSe (Ref. l32l ) 
with a small amount of Sn replacing Pb. A doublet is observed at 265 - 300 cm~^ which 
corresponds to structure in the sum spectrum of Fig. [S] (see (green) arrows) an is rather 
similar to the doublet described above for PbS (430 - 460 cm~^). The (green) arrows at 100, 
135, 175, 265 and 290 cm~^ also correspond to structure seen in the experimental spectra, 
although it is not clear whether the 135cm~^ structure is due to scattering by two phonons or 
forbidden scattering by one LO phonon.— A peak is also observed in the measured spectra 
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to difference 



at ~ 90 cm~^ (Ref. l3ll ). It has been attributed by the authors of Ref. 
scattering (LO-TO) although, according to FiglSl it could also be related to the peak in 
the sum DOS. Without measurements of the temperature dependence of this peak it is 
not possible to clarify the assignment. The Raman data available for PbTe seem to be of 
poor quality because of segregation at the surface of Te02 and other compounds: The only 



structure that is clearly identifiable corresponds to a single LO phonon (~120 cm ^).- 
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VI. HEAT CAPACITIES 



Figure displays the heat capacities of PbX (X = S, Se, Te) with the natural isotope 
composition on either constituent. For comparison, also literature data are shown. We have 
chosen to plot the measured Cp (and the calculated Cy ) so as to be able to read the data 
from the plots above ~ 50 K with some accuracy. Below ~ 50 K the heat capacities become 
rather small and it is more convenient to plot Cp^^/T^, as will be done in the subsequent 
figure. Notice that at the highest temperatures in Fig. [6], the heat capacities tend to the 
Petit-Dulong limit ( ~50 J/mol K). The agreement between experimental and calculated 
results is excellent although it is not possible, in these plots, to see the advantage of using 
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in the calculations a Hamiltonian with s-o interaction. Such advantage only appears below 
50 K, especially in the region where Cp/T^ has a maximum. 



60 




50 100 150 200 



FIG. 6: (color online) Temperature dependence of the heat capacities of PbX (X=S, Se, Te) 
from top to bottom, respectively measured by us and other authors. The results of our ah initio 
calculations, which include s-o interaction, are shown by (blue) solid lines. Literature data have 



been taken from Refs. 
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38, as indicated. 



Figure [7] displays the temperature dependence of C^^p/T'^ as measured in the 2 - 60 K 
range and calculated with and without s-o interaction for PbS and PbSe, all constituents 
with the natural isotope abundance. The data for PbS are identical to those published in 
Ref . y. 

It is of interest to compare the effect of the s-o interaction on the calculations of Fig. [7] 
with that found for bismuth (Ref. which has a s-o interaction similar to that of lead. The 
maximum in C^/T^ calculated for Bi without s-o interation is 20% below the measured one. 
Inclusion of s-o interaction yields a maximum 7% higher than the measured one, i.e., 27% 
higher than that calculated without s-o interaction.- The strong s-o effect calculated for Bi 
may result from the fact that there are two equal atoms per PC of Bi whereas for the lead 
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chacogenides there is only one heavy atom (Pb) per PC. 

Figure [7](b) displays Cp/T^ as measured for three samples of PbSe with different isotopic 
compositions, together with three curves for C^/T^ calculated with different contributions of 
the s-o interaction and with no s-o interaction. The purpose of this exercise was to identify 
the separate effects of the s-o interactions of Pb (atomic s-o splitting of the Qp electrons 
Apb= 1.27 eV) and Se (Ase =0.42 eV (Ref. Im))- The calculations in Fig. El^b) correspond 
to full s-o interaction, interaction only for Pb and interaction only for Se. Notice that the 
calculation with only the s-o coupling of Se is rather close to that with no s-o coupling at 
all. The s-o effect on is roughly proportional to the square of the s-o coupling parameter 
(cf. Ref. |6|) which for the valence electrons of Pb is three times bigger than that of Se. 
Hence, the effect of the s-o of Se alone is expected to be ~ 10 times smaller than that of 
Pb, i.e. unnoticeable in Figure [3 However, C„/T^ with s-o coupling for both constituents 
is slightly larger (~3%) than if only the s-o interaction of Pb is included. This suggests the 
presence of bilinear terms in the corresponding perturbation expression. In order to quantify 
these effects we define three partial contributions cpb , cse and cpb-so, which represent the 
contribution to the s-o effect quadratic in the s-o of Pb, in that of Te, and the bilinear 
contribution. A fit of these effects on Cu/T^ shown in Fig. [7|^b) yields: 

cpb = 99.6/iJ/molK^, csc ~ — S/iJ/molK'*, cpb-se = 76.2/iJ/molK^ (3) 

Obviously, the bilinear term cpb-so plays a rather important role in the perturbation 
expansion, almost as important as the quadratic term. This suggests the detailed inves- 
tigation of a material with a larger s-o coupling for the anion, such as PbTe (Atc =0.86 



eV, Ref. l34l ). The temperature dependence of Cp/T^ measured for samples with several 
isotopic compositions in the 3 - 14 K region, where the maximum occurs, is shown in Fig. 
IHl As already suspected, bilinear effects are large. They can be elucidated by multiplying 
the s-o Hamiltonian corresponding to the two constituents by two parameters Apb and Axe, 
respectively. Full s-o interaction is obtained for Apb = Atc= +1- Information about bilinear 
(and also cubic, cf. Ref. |6|) terms is readily obtained by reversing the sign of one of either 
Apb or Atc • Additional information concerning the perturbation expansion of Cp/T^ versus 
s-o interaction is obtained by setting one of either Apb or Axc equal to zero, the other equal 
to one. 

Figure [8] shows the experimental points for Cp/T^ of PbTe obtained with the isotopic 
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FIG. 7: (color online) (a) Heat capacity of PbS divided by T^. The (black) circles represent our 
experimental data (from Ref. ^). The (black) solid lines represents the ab initio calculations with s- 
o interaction for Pb and S, respectively (Apfe=A5=l). The red (dashed) line represents the ab initio 
calculation with s-o interaction left out {Xph=Xs=0). Note that for the calculations without s-o the 
maximum lies 27 % below the measured one. Inclusion of s-o interaction reduces the discrepancy 
by a factor of two. For the three curves (experimental and ab initio) the maxima occur at ~12.5 
K. (b) Measured Cp/T^ for three PbSe samples with different isotopic compositions as explained 
in the inset. The heat capacities of the three samples with different isotope composition are almost 
indistinguishable. The results of the ab initio calculations are represented by the (red) dashed, the 
(green) dotted and the (blue) and (black) solid lines, the latter two are almost indistinguishable, 
i. e. the effect of s-o interaction of only Se is not noticeable, while there is a significant increase 
near the maximum of Cp/T^ (Tmax ~ 11-3 K) if s-o interaction for Pb is taken into account ((red) 
dashed and green (dotted) lines. For more details see inset). 
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compositions given in the inset. The vertical scale is not wide enough to see differences for 
the various isotopes: such differences will be discussed later. The experimental maximum 
occurs at ~ 8.5 K, nearly independently of isotopic composition. In order to check the 
separate effects of Apb and Axc on the ab initio calculations of C.^,/T^ we have multiplied 
the s-o interaction of Pb and Te by Apb and Axc , respectively. We have then performed 
calculations for the six (after subtracting Apb = Atc = 0) sets of values of these separate 
parameters, as was done in the case of bismuth in Ref joj (only for a single A corresponding 
to bismuth, since there we dealt with a monatomic crystal). The seven sets of values of Apb 
and Atc used, also given in the inset of Fig. [HI are: (1,1), (1,-1); (-l)!); (-1)0) , (1,0), (-1,-1) 
and (0,0) for (Apb,ATe)- We have fitted the values of C^/T^ at 9K (the vertical line in Fig. 
[H]), after subtraction of the value for Apb = Atc) = (i.e., without s-o interaction for both 
constituents) to four parameters chosen as follows: 



a/T3(Apb, Atc) - a/T=^(Apb = 0, Atc = 0) = Xl^c, + X^^c^ + ApbAxeCs + A^bC4 (4) 

With five calculated values of Eq. (j4]) and four fitting parameters the corresponding values 
of Cj are overdetermined and only an approximate fit is expected. 

By minimizing the corresponding variance we find the values of Cj given in Table [I] to- 
gether with the calculated values of Cr^/T'^ (Apb ,ATe) and the fitted ones. The fit with four 
parameters of the five data seems to be reasonable and thus the perturbation expansion of 
Eq.(jl]) is acceptable. Like in the case of PbSe discussed above, the quadratic term ci is close 
to the bilinear one C3 . Also like in the case of Bi, with a s-o 6p splitting similar to that of 



Pb (Ref. l34l ). a cubic term is needed to represent the calculated results.- 



VII. DEPENDENCE OF THE HEAT CAPACITIES ON THE ISOTOPIC MASS 

The derivatives of C^JT'^ with respect to the isotopic masses of either Pb or S in PbS have 
been presented in Ref. |9|, obtained experimentally and by ah initio calculations. Although 
the calculations were performed without s-o coupling, the scatter of the experimental data, 
due to the small range of isotopic masses available, does not warrant a repetition of the 
calculations with s-o coupling. We shall present here similar results for PbSe and PbTe, with 
the ah initio calculations including s-o interaction. Figure [9] displays the results obtained 
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FIG. 8: (color online) (a) Heat capacity of PbTe divided by T^. The symbols represent our experi- 
mental data for various isotope compositions (cf. lower inset). The (black) solid line represents the 
ab initio calculations with s-o interaction for Pb and Te, respectively (Apb=ATc=l)- The (colored) 
solid lines represent the ab initio calculations with various combinations of the s-o interaction in- 
cluded (cf. upper inset, from top to bottom in the same order as in the figure). Literature data 



have been taken from Ref. 



381 ). as indicated. 



TABLE L Expansion, ciAp^ -|- C2A^g -|- c^XphXTe + (^i^pb used to quantify the contributions to 
the quantity Cp/T^(9K) if s-o coupling for either Pb (Ap6=l) or Te (AT'e=l) is taken or not taken 
(Ap6=ATe=0) into account (cf. Eq. The values of the coefficients Cj that optimally fit the data 
are: 01=247.5// J/molK'^, C2=64.2/x J/molK'^, C3=287.25// J/molK'^, and C4=-50.62// J/molK^. 



CsXpbXre + C4^pb 


ACp/T3(9K) (/i J/molK^) 
calcul. experim. 


Cl + C2 + C3 -1- C4 


548.3 


565.7 


Ci -1- C4 


196.9 


162.0 


C2 


64.2 


29.3 


Cl + C2 - C3 + C4 


-26.2 


-8.8 


Cl + C2 - C3 - C4 


75.0 


75.0 
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for PbSe using the isotopic samples given in the inset. 




FIG. 9: Derivatives of C^^p/T^ of PbSe with respect to the isotopic masses of (a) Pb and (b) Se 
as measured using samples of different isotopic compositions, compared with ah initio calculations 
for various isotope masses with s-o interaction included. 

A model for relating the position in the maxima of Figs. and [TU] was presented in 
Ref. y. It made use of a single Einstein oscillator to represent prominent peaks in the 
one-phonon DOS (Fig. |4]). The peak of the mass derivatives (Figs. [Qland fTOl) takes place at 
the Einstein frequency (converted into temperature by multiplying by 1.44) divided by 6.25. 
Using this procedure, we find for the maximum in Fig. [9]at 9 K {Mpb derivative) an Einstein 
frequency of 40 cm~^ which corresponds fairly well to the strongest peak seen in Fig. IHfor 
PbSe. For the Se derivative we find an Einstein frequency of 90 cm~^, in good agreement 
with the Se-like peak seen in Fig. |H Using the same procedure for Fig. [TU] we find from 
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the Mpb derivative the Einstein frequency of 28 cm~^, which also corresponds rather well to 
the sharp Pb-like peak of Fig. 4 for PbTe. For the Mpe derivative we estimate an Einstein 
frequency of 38 cm^^, which agrees with the maximum in the Te projection of the DOS in 
the corresponding vignette of Fig 4 (lowest band of PbTe). This contribution of Te to the 
mass derivative in the accoustic branch occurs because of the considerable contribution of 
the anion vibrations in this region, as displayed in Fig. HI A similar anion contribution is 
negligible in the cases of PbSe and PbS. 




FIG. 10: Derivatives of Cy/T^ of PbTe with respect to the isotopic masses of (a) Pb and (b) Te 
as measured using samples of different isotopic compositions, compared with ab initio calculations 
for various isotope masses with s-o interaction included. 

Since most of the Te weight to the DOS is in the optical phonon region (90 cm~^. Fig. 
m ) there should also be a corresponding peak or band in the lower vignette of Fig. [TUl 
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The experimental data in this region scatter too much to allow observation of this peak. 
The calculated curve, however, should exhibit some structure in this region. Although it is 
hard to see in Fig. [TU] it is possible to detect it if the derivative with respect to T of the 
mass derivative is plotted. This is done in Fig. [TTl In the lower vignette of this Figure, 
concerning the M^e , we see indeed weak structure between 20 and 40 K which corresponds 
to the upper PbTe band of Fig. 1 (gOcm'^x 1.44/6.2= 21K). 




FIG. 11: (color online) Derivative with respect to T of the calculated curve of Fig. [10] (lower 
vignette), showing a weak band in the 20-40 K region which corresponds to the Te-like band in 
Fig. II for PbTe. See text. 



VIII. TEMPERATURE DEPENDENCE AND ISOTOPIC MASS DEPENDENCE 

OF Cv,p/T^ 

As already mentioned, the dependence of Cy/T^ and Cp/T^ on isotopic mass, or more 
precisely the corresponding derivatives, can be obtained from the temperature dependence 
of Cy^p/T^ in the case of monatomic materials. The corresponding relation for diatomic 
materials has been given in Eq. ([2]). The correctness of this equation has been demonstrated 
in Fig. 10 of Ref. l9| for PbS. In Fig. [12] we present one further example corresponding to 
PbSe. 

The (red) dashed-dotted curve represents the derivative of the calculated logarithmic 
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FIG. 12: (color online) Illustration of the relationship between the temperature dependence of 
Cv,p/T'^ and its dependence on the isotopic masses of both constituent atoms, drawn from ah 
initio calculated data for PbSe with natural isotope abundance for Pb and Se, respectively. The 
(black) circles represent experimental data for the temperature dependence {dln{Cy^p/T^) /dlriT), 
the (black) solid line represent the l.h.s of Eq. ([2]). The other quantities are explained in the inset. 

derivative of C^/T^ versus T. The green and blue curves represent the corresponding isotopic 
mass derivatives as found in Eq. ([2]) (obtained from the calculated data (Apb = Atc = !)• 
The experimental data scatter too much to obtain meaningful results). The black curve is 
the sum of these derivatives. According to Eq. ([2]), the black and red curves should be 
indistinguishable. This is indeed the case for T > 30K. At lower temperatures the scatter 
of calculated and measured points becomes rather large because both C^^p and are very 
small. Nevertheless reasonable agreement is found. It is interesting to note that a change 
in the Pb mass strongly affects the T dependence of C^JT^ around 5 K whereas the change 
in the Se mass affects that T dependence in the 30-50 K region. 

IX. CONCLUSIONS 

We have presented theoretical and experimental investigations of some lattice properties 
of three lead chalcogenides (PbS, PbSe, PbTe) with rock salt struture. Because of the 
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presence of lead, these materials provide an excellent case study for the elucidation of the 
effects of spin-orbit interaction on lattice properties as derived from ab initio electronic 
structure calculations. Among other properties that have been investigated are the lattice 
parameters, the phonon dispersion relations, the heat capacity and their dependence on 
isotopic masses. The effect of spin-orbit interaction on the lattice parameters is less than 
1% and can be neglected compared with other sources of computational error. Its effect on 
the lattice dynamics and the corresponding low temperature heat capacities is considerable. 
For these properties spin-orbit interaction significantly improves the agreement of theory 
with experiment. Because of the large scatter in the effect of isotope masses on the heat 
capacity, related to the small range of of isotope mass available, spin-orbit interaction does 
not play an important role in the comparison of theoretical and experimental data for the 
isotope mass dependence of the heat capacity. This dependence can be modeled on the basis 
of two Einstein oscillators. 

We have calculated, from the ab initio phonon dispersion relations, the corresponding 
density of one-phonon states and those of two-phonon states relevant for optical spectro- 
scopies. The latter have been compared with few extant second order Raman spectra of 
these materials. 
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